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ABSTRACT 

A major uncertainty in accretion disk theory is the nature and properties of gas turbulence, which drives transport in protoplanetary 
disks. The commonly used viscous prescription for the Maxwell-Reynolds stress tensor gives rise to a meridional circulation where 
flow is outward near the midplane and inward away from it. This meridional circulation has been proposed as an explanation for 
the presence of high- temperature minerals (believed to be of inner solar system provenance) in comets. However, it has not been 
observed in simulations of magnetohydrodynamical (MHD) turbulence so far. In this study, we evaluate the extent to which the net 
transport of solids can be diagnostic of the existence of meridional circulation. To that end, we propose and motivate a prescription 
for MHD turbulence which has the same free parameters as the viscous one. We compare the elfects of both prescriptions on the 
radial transport of small solid particles and find that their net, vertically integrated radial flux is actually quite insensitive to the flow 
structure for a given vertical average of the turbulence parameter a, which we explain. Given current uncertainties on disk turbulence, 
one-dimensional models are thus most appropriate to investigate radial transport of solids. A corollary is that the presence of high- 
temperature material in comets cannot be considered an unequivocal diagnostic of meridional circulation. In fact, we argue that 
outward transport in viscous disk models is more properly attributed to turbulent diffusion rather than to the mean flows of the gas. 

Key words, accretion, accretion disks ~ instabilities - turbulence - magnetohydrodynamics (MHD) ~ meteorites, meteors, meteoroids 
- comets 



1. Introduction 

Extraterrestrial samples provide ample evidence for signifi- 
cant radial transport in the protoplanetary disk from which the 
solar system emerged. For instance, carbonaceous chondrites 
contain high-temperature materials such as refractory inclusions, 
believed to have formed close to the Sun in the first stages of 
disk building (|Wood,,2004^ , chondrules, formed 1-3 Ma later 
( [Villeneuve et al.||2009[ l, as well as aqueous alteration products 
indicating the presence of ice when these meteorites accreted. 
Yet more spectacularly, rare refractory objects and chondrule- 
like fragments have been identified in dust returned from comet 
Wild 2 ( jZolensky et al.|2006,, Bridges et al.|20T2| . 

One way of transporting solids outward from the vicinity of 
the Sun is through bipolar outflows, as in the X-wind model of 
I Shu et aL] ( [2001 j (see also Hu„2010j ; however, whether solids 
can be present at and efficiently transported from the X point has 
been called into question by Desch et al. (2010 1. More recently, 
scenarios have been proposed where transport is effected inside 
the disk by turbulent motions of the gas — that same turbulence 
which would also drive accretion of the disk gas onto the central 
star ( Bockelee-Morvan et al.|2002 



2005';'Johansen & Kla hr|2005[|Cies 



Boss 2004; Cai-balUdo et al. 



2009 i Hughes & Armitage 



2010; Jacquet et al. 2011 Yang & Ciesla 2012). hi early stages. 



where the disk is massive, such turbulence could be due to gravi- 
tational instabilities ( Boss|2004 1 whereas magnetohydrodynam- 



ical (MHD) turbulence powered by the magnetorotational insta- 
bility (MRI; Balb us & Hawley|1998) l should dominate afterward. 
However, the MRI is expected to be suppressed over a significant 
range of heliocentric distances because of insufficient ionization 
( |Gaiimiie]|1996| l, so that in this region, referred to as the "dead 



zone", turbulence should be reduced and have a more hydrody- 
namical character 

While many studies have restricted attention to one- 
dimensional disk models, where only the vertically averaged tur- 
bulence parameter a had to be prescribed, some have attempted 
to calculate the transport of solids in two- or three-dimensional 



disk models (e.g. Hersant et al. 2001] Ciesla 2009 ; Charnoz et al 



201 1 Ciesla|2010a 201 1 1. Going to 2D or 3D requires however 



more assumptions on the properties of turbulence to be made de 
spite the considerable uncertainties of accretion disk theory and 
the diversity of potential turbulence drivers alluded to in the pre- 
vious paragraph. A widespread prescription models turbulence 
as an effective viscosity in the averaged dynamical equations 
of the gas (e.g. Takeuchi & Lin| 2002): the circulation resulting 
from this "viscous" prescription typically involves an outward 
flow around the disk midplane and inward flows in the "atmo- 
sphere", with the vertically integrated flow being inward (i.e. a 
positive mass accretion rate). This meridional circulation was 
found by Ciesla ( 2007) to improve outward transport of inner so- 
lar system material and hence account for the presence of high- 
temperature minerals in comets. However as yet, the meridional 
circulation is an essentially theoretical construct that has not 
yet been observed in numerical simulations of turbulent disks, 
in particular in the global MHD simulations of Fromang et al. 
( 201 1\ and [Flock et al. ( 2011[ l. As yet, the vertical structure of 
the gas flow is uncertain. 

The purpose of this article is to investigate, through analyti- 
cal calculations, to what extent the net outward transport of small 
solids can be considered a diagnostic for meridional circulation. 
We propose and motivate a prescription for the vertical profile 
of the Maxwell-Reynolds stress tensor in the case of MHD tur- 
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bulence and calculate the resulting mean flow of the gas. The 
prescription, inspired by previous numerical studies, formally 
depends on the same parameters as the standard "viscous" pre- 
scription, allowing a direct comparison. In particular, we com- 
pute the net (vertically integrated) radial flux of solids in both 
prescriptions. It is found that this net flux is weakly dependent on 
the prescription used, which will be explained, and that therefore 
meteoritic and cometary properties cannot be used as evidence 
for a specific vertical flow profile. 

This article is organized as follows: in Section |2] we intro- 
duce general notions on gas turbulence and the two prescriptions 
considered here. In Section [3] we investigate and compare the 
dynamics of solids in both prescriptions. We discuss the results 
in Section |4]before concluding in Section |5] 

2. Modeling of gas turbulence 

In this section, we introduce the modeling of gas turbulence. 
After having reviewed general equations, by way of establish- 
ing notation, we introduce the two prescriptions of the Maxwell- 
Reynolds tensor investigated in this paper, namely the viscous 
and the MHD prescriptions. 

2.1. Generalities 

The disk is described in a cylindrical coordinate system, with 
R the heliocentric distance, z the altitude above the midplane, 
and the azimuthal angle. We note u, p, T and P = pc^ the gas 
velocity, density, temperature and pressure, respectively, with c, 
the isothermal sound speed, = -\/GMq/R and Q. = v^/R are 
the Keplerian linear and angular velocities, respectively. 

We treat the disk as vertically isothermaQso that vertical hy- 
drostatic equilibrium implies the following density stratification: 



^/2^H(R) 



exp - 



2H{Ry 



(1) 



with the pressure scale height H = cJQ. and the surface density 

We assume the disk to be turbulent and assume axisymmetry 
in the sense that variations of any quantity Q in the azimuthal 
direction may be treated as turbulent fluctuations about a mean 
Q. We denote the Eulerian perturbations with 6Q = Q - Q. For 
any quantity Q and a weight function w, we also define the w- 
weighted vertical average of Q as: 



<G>vv = 



rr Qiz)w(z)dz 



(2) 



Because of the pressure gradient and the departure from 
the equatorial plane for z + is not exactly equal to the 

Keplerian velocity vk- Indeed, centrifugal balance yields: 



1 dP 



'ii+izm'" 
1 



pdlnR 



^ Vk + 



dP(R,0) ndlnT ^ llH^^ 
z + 



2p(R,0)n dR 4 dR 



' Although this approximation breaks down in the surface layers of 
the disk (e.g. [Chiang & Goldreich|1997 1, this is of no concern here as 



this generally affects a small fraction of the surface density of the gas 
and an even smaller one of that of the (partly settled) solids. 



where we have used z,H<kR. Since, generally, dP(R, 0)/dR < 
and dlnT/dR < 0, the flow is subkeplerian (e.g. Fromang et al. 
Flock et al.|201 1| . However, the approximation 



EOTT 



Vk 



will be generally sufficient except where derivatives in z will be 
invoked (in Section [2!2] l. 

The angular momentum equation, averaged over turbulent 
fluctuations, reads ( [Balbus & Papaloizou|1999] ): 



dRu0 1 d 

puR— 1 — — 

^ dR RdR 



[r^Tr 



(4) 



where Tr^ and T,^, the R(p and z0 components of the turbulent 
stress tensor, are defined as: 



Ti^ = p6u^5ui 



(5) 



with B the magnetic field and / = R, z- We parameterize Tr^ as: 



Tr^(R,z)= -a(R,z)P(R,z). 



(6) 



Note that equation (j6]) is only a definition and does not presup- 
pose any prescription of the stress tensor. In particular, a may a 
priori vary with R and z. The arbitrary factor 3/2 has been intro- 
duced to facilitate direct comparison with the turbulent viscosity 
formalism. 

Equation d?]) may be rewritten as: 



Ur = -- 



1 d IR"'- 



R^'^dRX Q 



1 dT,. 



n dz 



(7) 



From now on, we will drop the overbars. The latter equation may 
be in tegrated to yield the vertical d ensity- weighted average (mr)p 
(e.g. |Lynden-Bell & PringlellOTil P 



{ur). 



PUrAz - - 



m^i^ dR 

3 d 

'm^i^dR 



R^'^I.{a)p 



TR<i,dzj 

11. 



(8) 



Thus, the net radial transport of gas depends only on {a}p as far 
as turbulence properties are concerned. 

We now turn to the two prescriptions for the Maxwell- 
Reynolds tensor to be considered in this paper, namely the "vis- 
cous" and the "MHD" prescription. 

2.2. Tine viscous prescription 

The viscous prescription, first introduced for hydrodynami- 
cal turbulence by Boussinesq ( 1 877 1 and Reynolds ( 1895) — but 
see e.g. Schmitt ( 2007|l — ~and co mmonly used in protoplane 



tary disk contexts (e.g.'Urpin 1984 Takeuchi & Lin 2002 |Ciesla 
2009 ) , assumes that turbulence can be modeled as an eiffective 
kinematic viscosity v, so that: 



W.) = -pv4(|) 



pM^-pv—^-pv-^ 



Q51nr 



dR 



(9) 



(10) 



where we have made use of equation 

' This assumes that T-^ vanishes at infinity, which may not hold in 
case of steady outflows, but even then, its contribution to {ui{)p is likely 
negligible jBai & Stone|2012 Fromang et al. 2012 1. 
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Gas velocity profile 




numerical studies ( 


Miller & Stone|( 


2000|l; 


Hirose et al.|(|2006|l; 


Flaig et al. 
et al. (2010 


(•2010 


; Fromang & ^ 


elson ( 


2006|l; 


Dzyurkevich 


;,Sorath 


ia etal. (2010|l; 


Promang et al. ( 


2011 1; Flock 


et al.|(|20li 


l;|Guan & Gamnue|(|2011|l;|Fromang et al. (|2012|i; 


Bai & Stone (2012|l; see also Uzdensky (2012|l), the vertical pro- 



file of Tr^ shows a plateau for |z| lower than a given z^^^. It then 
falls off, propoitional to the density, in the corona, where mag- 
netic pressure is comparable to the gas pressure, quenching the 



MRI and leaving a transsonic turbulence there (e.g. Flock et al. 
|2011j). We thus adopt: 



- z 



2//2 



e(\z\ 



Fig. 1. Plot of the mean radial velocity of the gas for the viscous 
(dashed) and the MHD (continuous line) prescriptions. We have 
taken E, T oc R^^ ''^, and a{R,Q) oc R" (with the local value of 
{a)p being 10"^), with three values of the exponent a used in the 
plots. The viscous prescription entails a parabolic flow profile 
(the so-called meridional circulation) while the MHD prescrip- 
tion entails an inverted gaussian profile (until the corona, here at 
kl > Zmax = 3.0//). Here, for a = -1, the vertically integrated 
flow is outward in both prescriptions. 



Comparing equation (j9]l with equation (j6]) leads to the iden- 
tification: 



rR//?,z) = r«^(/;,0)exp( 

with 6 the Heaviside function, defined as: 

..,^_/Oif-^<0 
^W-|iif_^>0 

We thus have 

a(/?,0)e^'/2H^if|z| <Z„,x 



a{R,z) = 
with 



if kl > Zn 



= //, 21n 



a{R, 0) 



(14) 



(15) 



(16) 



(17) 



As to ffmax, we will adopt in numeiical applications the fiducial 
value 1/3 tha t aiises in the ideal MHD limit of the study of Bai &| 



Stone (201 1|| This clearly differs from the viscous prescription. 



v{R,z) = a(R,z) 



(11) 



We hence have: 



{a)p^ —a{R,Q). 



+a„ 



vCrfc 



(18) 



At this point the viscous presciiption is non-tautological only in 
the presciiption of T,^. It is generally further assumed that a is 
constant, at least vertically, i.e.: 

a(R,z)^a(R,0)^{a}p. (12) 
The radial velocity may then be calculated from equation (|7]l: 



ur(R,z) 



a(R, 0)c2 



Q 



d_ 
dR 



i„((«(«,o,a'H).A,„g)(£) 



with erfc the complementary Gauss error function defined as: 

2 r°° . 

erfc(.x) = 1 -erf(x) = — e'^'dy (19) 
yn Jx 

Note that we have made no assumption on the value of a(R, 0) 
or on how it scales with other disk parameters as we focus here 
on the z dependence of the flow. 

Considerably less is known about T.^. We propose a break at 
(li^l ~ -^max similar to that of Tr^ in the form: 



This flow, parabolic in z, is known as the meridional circulation 
(see Fig.[T]l. Typically, the velocity near the midplane is positive 
( Takeuchi & Lin|2002[ l, but the inward flows away from this re- 
gion result in a positive net (vertically integrated) mass accretion 
rate (i.e. toward the Sun; see equation (|8]l). This circulation may 
be interpreted as follows: at the midplane, the large radial den- 
sity gradient (p oc 1,0./ Cs) makes the viscous torque exerted from 
the inner, more rapidly rotating disk regions greater than that 
received from the outer, slower rotating disk, hence a net gain 
of angular momentum and an outward flow; at high altitudes, 
however, the density gradient is reduced and even switches sign, 
changing the direction of the effect (e.g. Takeuchi & Lin|2()02 l. 



T,4R,z) = C(R,z)exp 



2//2 



-0(k|-Zmax) , 



(20) 



2.3. MRI-turbulent disk prescription 

Simulations of magnetohydrodynamic turbulence have so far 
failed to exhibit a meridional circulation |Fromang et al.|2011[ 
Flock et al. 201 1) and this may be traced to the discrepancy be- 
tween the measured stress tensor and that assumed in the vis- 
cous prescription ( [Fromang et al.|20TT] l. As noted by numerous 



with C(R,z) an as yet unspecified smooth function which sym- 
metry about the midplane requires to be odd in z- By requiring 
that ur as expressed by equation (j7]i be continuous at z = +Zmax, 
we have: 

C(R, z„ax) = -C(R, -z.™x) = Tr^(R, 0) ^ , (21) 

which essentially sets the scale of T^^ ~ (////?)/'«^|^If we adopt 
the simplest dependence C(R, z) oc z, we obtain: 

T.^,(R,z)^Tr,(R,z)^^z (22) 

^ We have recast their equations 23 and 25 in terms of our definition 
of a (equation (j6|); a^..,^ corresponds here to the limit where gas and 
magnetic pressure are equal. 

Note that while our derivation makes T-^ depend on the formally 
nonlocal quantity dzmnxIdR, this is not inconsistent with the fact that lo- 
cal simulations can exhibit a nonzero T^^ (see e.g. [Fromang et al.|2012[ l 
to which the radial gradient of a might adjust. 
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This yields a profile similar to that measured by |Fromang et al.| 
( |20n] l. 

Equations 



14 1 and ( 22 1 constitute what we will refer to as 



the "MHD prescription". It is noteworthy that it involves the 
same free parameters (a(R, 0), or equivalently {a}p, and its pos- 
sible radial derivative) as the viscous one. 

The mean radial velocity resulting from this prescription 
may then be calculated. For |z| < Zmax, 



ur(R,z) 



d 



dR 



1/2 



Tr4R,Q) 



Mfi(0)e37^,(23) 



and for |z| > z„ 



ur(R,z) = 



d IR^'^ 



d_ 
"dR 



1 / ^max \ / ^ \^ 



(24) 



The flow is plotted in Fig.[T] In contrast to the meridional circula- 
tion, velocities do not change sign (except perhaps in the corona 
if da(R,0)/dR + 0). Indeed the stress tensor no longer scales 
with the density vertically, in contrast to the viscous prescrip- 
tion, and thus the effect of the change in radial density gradient 



with height (see Section 2.2 1 disappears. Then, the radial veloc- 



ity does not switch sign, and always has that of its average (at 
least outside the corona). The profiles are similar to those aris- 
ing from the global simulations of Fromang et al. (201 1 1. They 



differ, however, from those of |Flock et al.| ( |201 l| l, which resem- 
ble an inverted meridional circulation profile, where the posi- 
tive velocities in the corona correspond to outflows potentially 
escaping the disk. Such outflows launched from the disk have 
been also recently studied in local simulations with uniform ini- 
tial vertical magnetic field ( |Lesur et al . 2012; Bai & Stone 2012; 
Fromang et al. 2012 1 and if present, are thus not captured by 
our prescription of T^^ (which would not vanish at infinity for 
outflow-launching simulations). We shall however argue later 
that its exact form has little effect on our conclusions. Our pre- 
scription for the Maxwell-Reynolds tensor is thus only intended 
here to provide a sensible MHD counterpart to the viscous one 
for comparison purposes. Certainly, more measurements of Tr^ 
and T,fj, in stratified numerical simulations are needed to arrive 
at a more definitive prescription. 

3. Impact on the dynamics of solids 

With the two above prescriptions on the gas flow in hand, we 
now investigate their impact on the dynamics of solids tightly 
coupled to the gas. We first briefly review the basic concepts of 
aerodynamic transport of solids in disks, before focusing on the 
vertical distribution of particles and their net radial motions. 



3.1. Generalities 

Let us consider a population of identical spherical particles 
of internal density and radius a embedded in the gas. We de- 
note by pp the density of this population (viewed as a fluid) and 



Ppdz 

CO 



(25) 



their vertically integrated surface density. We also define a nor- 
malized solid-to-gas ratio f by: 

PpiR,z) _ ^p(R) 



p{R,z) m) 



fiR,z). 



(26) 



The effects of gas drag is characterized by a stopping time 
T (see e.g. Weidenschilling|[T977[ ) from which one can define a 
measure of the coupling of the particles to the gas on an orbital 
timescale as 



St = Q.T, 



(27) 



which we will take to be <s; 1 in accordance with our tight cou- 
pling assumption. 

Neglecting any feedback of the solids on the gas, the conti- 
nuity equation averaged over turbulent fluctuations reads: 



dpp 



1 d 



dt ^ RdR 



d_ 
'Jz 



R \PpVp,R - DrrP 

PpVp,z - D,,p 



d_IPp 
dR\p 
d_lpp 
dz \ p 



(28) 



where we have introduced the mean velocity Vp of the particles 
given by ( Youdin & Goodman|2005| l 



V,, = U H- Vdrift = U + T- 



(29) 



and also diffusion coefficients Drr and D , parameterized as: 



D, 



rr 



■Q 



and D,- 



Q' 



(30) 



with 6r and 6, dimensionless parameters of order a. We will 
assume that the radial and vertical Schmidt numbers 



Scr= — 
Or 



and 



(31) 



are vertically constant. 

Following [Jacquet et al. ( 2012[ l, we coin: 



S = 



St 



^drift.R 
Ur 



and S, = 



St 



5 X Sc- 



(32) 



which measure the coupling between particles and gas on a 
global scale, in the radial and vertical directions, respectively. 

3.2. Vertical distribution 



Equation ( [28) is dominated by the vertical flux balance be- 
tween settling and diffusion. On timescales longer than the ver- 
tical mixing timescale 



1 



Qmax((J,,St)' 



the vertical distribution of the particles obeys: 
^--S-SR,z)-,, 



(33) 



(34) 



where the role of 5, as a measure of settling is apparent. It may 
be noted that our assumption of constant Scj makes the vertical 
distribution independent of the prescription of T,^. 

In the viscous prescription, where a, and hence 6, is verti- 
cally constant, 5j oc p"' so that equation (34 1 may be integrated 
as ( jTakeuchi & Lin|2002| i: 



fiR,z) = 



ItiH 



£ J exp {-S,{R, z') - z'V2H^)dz' 



exp{-S,(R,z)) 



exp(-S,{R,Q)e2H2), 



(35) 
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Fig. 2. Plot of the dust-to-gas ratio normalized to the column 
density ratio (/) as a function of z, for both viscous (dashed line) 
and MRI-turbulent prescriptions (solid line) and for three values 
of the "settling parameter" S~ = St/((J;)p, with {S;)p - \Qr^. In 
the viscous prescription, this ratio falls off as a "double gaussian" 
whereas it only does so as a gaussian for the MHD prescription 
over the bulk of the disk thickness, before decreasing as a "dou- 
ble gaussian" in the corona (here at |z| > z^ax - 3.QIT). 

a "double gaussian" which drops rapidly at large \z\. 

In the MHD prescription, S~ is vertically constant for |z| < 
Zm-dx before increasing as p"' in the corona. This gives: 



fiR,z)=mQ)sxpl-S,{R,Q) 



for |z| < Zmax, and 

VS;(«,0) 



f(R,0) 



laiR,0) 



exp 5;(;;,0) 1 - 



a(R,Q) 



\ ^max 

for \z\ > Zmax, with 

erf( V(Sj(i?,0)+l)ln(ff„ax/Q'(«,0))) 



(36) 



(37) 



2_ ( a(R.O) \ 



y 5.(^,0)+! 

/^«0«-p(-/ + 5,(/^,0)(l - ^.-))dy|38) 

At first the solid-to-gas ratio falls off only as a gaussian (because 
of increasing a counteracting looser coupling due to decreasing 
p) and only in the corona does it decrease as a "double gaus- 
sian" similarly to the viscous prescription case. Thus, a purely 
gaussian fit tends to overestimate, and a "double gaussian" one 
to underestimate the dust density away from the midplane in this 
case, as observed by Fromang & Nelson (2009) in their simula- 
tions of dust settling in MRI-driven turbulence. 

The vertical distribution in both prescriptions is plotted in 
Fig. |2] Similarly to |Fromang & Nelson (20091, we suggest that 
observations of present-day protoplanetary disks, in probing dust 
remaining on their surface, could discriminate between the two 
vertical distributions. 



3.3. Net radial flow 

We now turn to the radial motions of the particles 
Specifically, we want to integrate the continuity equation ([28] 



Fig. 3. Plot of the net vertically averaged velocity Vp id as a func- 
tion of S - St(R,0)/{a}p for the viscous prescription (dashed 
line) and the MHD prescription (solid line), for three values of 
Sc, (0.1 in blue, 1 in green and 10 in red). The velocities are 
always inward and we plot here their absolute value. We have 
taken E, T cc R-°^\ da{R,0)ldR ^ Q, 6r = 6, and {a}p = 10"^. 
The dotted line corresponds to the simple one-dimensional ap- 
proximation (see equation (42 1) which is common to both pre- 



scriptions and independent of Sc^. The approximation is arbi- 
trarily good in the limits 5 <K 1 and 5 » 1 and the imparted 
error does not exceed a factor of two when 5^1 except for the 
viscous prescription at Sc-.=10 where it can reach one order of 
magnitude. 

over z and deduce the net radial flow for both prescriptions. 
Formally, the result of this integration reads: 



5Ip 
dt 

with 



RdR 



dR \ Z 



0, 



where 

dfpAz 



(39) 



(40) 



(41) 



is a correction from the diffusion term, which is generally nega- 
tive (see e.g. equation ( |A.4| i), as the disk thickens with increasing 
heliocentric distance, inducing its "surface" to diffuse particles 
inward into the rarefied, particle-depleted gas. 

If we restrict interest to timescales longer than fvm, we can 
adopt the equilibrium vertical distributions of the previous sub- 
section to calculate the above. Results are presented in appendix 
[A|and Vp id is plotted in Fig.[3]for different values of the vertical 
^hmidt number. The main conclusion to draw from this plot is 
this net velocity is weakly dependent on the exact flow structure 
and is generally well approximated by 

Vp,iD(-R) « {ur)p + VinmiR, 0), (42) 
a widespread formula in one-dimensional calculations (e.g. 



Cuzzi et al.|20()3||Ciesla|2070 F:'Hughes & Armitage' 20T0l[^ng| 
& Ciesla|2012 1. Agreement is arbitrarily good when S is either 
very large or very small, and the deviation is largest for 5^1. 
The deviation is generally less than a factor of 2, except for the 
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viscous prescription // the Schmidt number is large. But even 
then, the net radial velocity does not (in general) become posi- 
tive, as had been noted by T akeuchi & Lin, (2002) , even if they 
did not include the diffusion correction. 

Why does the ID approximation (|42|) provide a so good 



Inward transport zone 



match, at least when Sc^ ~ Sc^ ~ 1? This stems from the facts 
that (see also [Jacquet et al.|2012| l: 

(i) In the limit S,S~ <K 1, particles are well-mixed with the 
gas and only drift negligibly with respect to it, so that their 
vertically-averaged velocity is indistinguishable from that of 
the gas, {uR)p. 

(ii) In the limit S,S^ » 1, particles are concentrated at the mid- 
plane and their motion is dominated by drag-related drift 
(evaluated in this region), so that their net radial velocity is 
essentially Vdi-ift(^, 0). 



As the expression (42i captures both limiting behaviors- 



with diffusion correction Vpcon being negligible in both limits 
(| Jacquet et al.|201^ — , it indeed provides a suitable approxima- 
tion for all possible values of 5 . It is noteworthy that the above 
reasoning is general and extends to prescriptions other than the 
ones considered in this study, e.g. if one adopts an alternative 
mathematical expression for T-^ in the MHD prescriptionPl 
In detail, some deviation is expected (and seen in Fig.ji]) for 

5 ~ I which is outside either above limit, but as the approxi- 
mations just begin to break down there, the order of magnitude 
at least should be accurate. If however Scj. is very different from 
unity, the domain of non-validity of the above regimes is wider, 
and depending on the vertical variations of the velocities, this 
could impart a more significant error. This is what we see for the 
viscous prescription if Sc, » 1. Indeed, under such conditions, 
one can have S < 1 < S^: then, the particles are concentrated 
around the midplane (since S , > 1) so that their gas velocity 
has an important contribution from the outward-directed flows, 
which is incompletely compensated by the inward drag-related 
drift (since S ~ Vdrift/fs < 1). hence a significantly less negative 
Vp.iD- We note however that large values of the Schmidt number 
are not expected in the hydrodynamical turbulence which the 
viscous prescription is intended to model ( ^Prinn, 1 990^ |Dubrulle 

6 Frisc h|1991 1, in contrast to MHD turbulence ( CarbaUido et al. 
'2005nJoha'nsen & KIahr|2005HJohansen et al.|2006| ), so that this 



situation appears unlikely anyway. 



4. Discussion 

From the preceding section, it appears that the net radial 
transport of particles is weakly sensitive to the vertical flow 
structure and is well approximated by the one-dimensional ex- 
pression given by equation (42 1. In particular, ;/ the mass accre- 
tion rate is positive (toward the Sun), the non-diffusive contribu- 
tion to the net radial flux should be generally inward. 

This might seem paradoxical in the case of the viscous pre- 
scription, as then, a population of small (S 1) particles around 
the midplane should be transported outward, even in the absence 
of turbulent diffusion. The paradox lies in our calculation per- 
taining to timescales longer than the vertical mixing timescale 
fvm (to warrant use of the equilibrium vertical distribution of 



The contribution to {ur)p of T,^ in our prescription is: 



c] d\nz„ 



n dR 

zero for 5 <c 1 and -2)(a{R, Q)c]lQ.)d\nzm:mldR <c <vj,-ift>p„ for 5 » 1. 



Vertical mixing 



e 



n-1 n n+1 

Outward transport zone 

Fig. 4. Cartoon of the toy model for meridional circulation. We 
distinguish between an "outward transport zone" (rightward ar- 
rows) and an "inward transport zone" (leftward arrows), which 
are assigned fractions 9 and \ - {9 < 1 /2) of the local col- 
umn density, respectively. Integers n-\,n,n+ \ represent radial 
locations whose separation corresponds to one vertical mixing 
timescale. Vertical mixing is symbolized by double vertical ar- 
rows. 



particles) and thus to radial scales larger than the corresponding 
radial excursion 

(H \ 
^===,— j<//.(43) 

So while a finger of material lying near the midplane would 
undergo some outward transport in the viscous prescription, it 
would not travel significantly further than this we ignore ra- 
dial diffusion. This we can see perhaps more concretely with a 
simple toy model, which we now present (see also Fig.|4]l: 

We consider the disk thickness to consist in an "outward 
transport zone" (the midplane), comprising a fraction 6 of ma- 
terial, and an "inward transport zone" (the "atmosphere"). We 
ignore radial diffusion and consider that both zones transport 
material at equal and opposite velocities, so that the net flow is 
inward if 9 < 1/2, as we shall henceforth assume. We discretize 
time and space (in the radial direction), with the temporal and ra- 
dial steps being fvm and AR^^^, respectively. Each time step cor- 
responds to two successive motions: advection one radial step 
forward or backward within each zone, and vertical mixing be- 
tween the two zones. Thus, if we denote by u„ t the total column 
density of a contaminant at time t and at radial location n — ^both 
natural integers, in units of fym and ARy^ — , Mn,t+i is given by: 



^n,t+l - ^^Mn-l,t + (1 ^ S)U 



n+l,t- 



(44) 



As boundary condition, we set uq , - 1 (in some arbitrary unit) 
and as initial condition m„ n = for n > 1. 



One can demonstrata°l that M„ t determined by equation (44i 



and its boundary and initial conditions satisfies: 

(i) Mn,t is a monotonically decreasing function of n. 

(ii) Mn,t is a monotonically increasing function of f. 

(iii) < u„ i < (0/(1 - 9))" (with the upper bound being a station- 



ary solution of equation (44 1) 



(iv) Unt converges toward {9/(1 - 9))" as f ^ +oo. 

Since 9 < 1/2, (0/(1 - 9))" and a fortiori Mn.t vanish as 
n — > -i-oo and in practice become negligibly small after several 
radial steps, corresponding, if we revert to physical units, to a 
radial excursion much smaller than R. The conclusion from this 
toy model is thus that no significant outward transport can be 



^ By recurrence on t for (i)-(iii) and noting that (ii)-(iii) implies con- 
vergence of u„ i and that the limit has to be a stationary solution of the re- 
currence equation, but must have zero projection on the solution (l)„eM- 
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achieved without including radial diffusion, even in a meridional 
circulation context. 

However, 2D disk models assuming the viscous prescription 
did include radial diffusion, and did achieve significant outward 
transport (e.g. |Ciesla|2009| l, but from the above, the latter should 
be viewed as being caused by the turbulent diffusion rather than 
the meridional circulation itself. Nonetheless, it does hold that 
the meridional circulation, in reducing somewhat the inward ad- 
vective velocities (see Fig. |3]l which counteract turbulent diffu- 
sion, enhances outward transport relative to the one-dimensional 
approximation, and this can have nonnegligible effects for long- 
range tran sport. F or example, in the "outward transport" sim- 



ulation of 



Ciesla 



(,2009) with M = 10"^ 
is obtained (after 10^ 



Mo/yr, a crystalline 
yr) at 10 AU, where 



fraction of 40 

the one-dimensional approximation would have predicted 17 % 
(using equation 29 of Jacquet et al. ( 2012| l). However, the one- 



dimensional model would have retrieved a crystallinity of 40 % 
had we taken Scs - 0.5 instead of the Sc« - 1 chosen by Ciesla 
( 2009| l, which is at least equally plausible. Thus, the differ- 
ences between meridional circulation and the one-dimensional 
approximation circulation are actually within the errors of tur- 
bulence parameters. A corollary is that the occurrence of high- 
temperature materials in the samples returned from comet Wild 
2 ( Zolensky et al.|[2006 1 is not diagnostic per se of meridional 
circulation. More generally, net radial transport of early so- 
lar system material is no compelling constraint on the verti- 
cal variations of the gas flow. The question of the form of the 
Maxwell-Reynolds tensor in protoplanetary disks, be it in MRI- 
active regions or in the dead zone, is thus most appropriately ad- 
dressed by numerical simulations of gas turbulence, although, as 
we mentioned earlier, observations of protoplanetary disks could 
also offer important diagnostics linked to settling. 

We caution before closing this section that we have not con- 
sidered, in this study, the possibility of outflows launched from 
the disk (e g. (Flock et al.||2011 ; Bai & Stone 2012; Fromang 
et al.||20l"2 Lesur et al.||2012| l~which would yield a sink term 
on the right-hand- side of equation ( 39 1 — potentially important 



for vertically well-mixed particles {S. <K 1). The possibility and 
quantitative importance of such outflows has yet to be explored 
in more details before their impact on the radial motion of solids 
can be satisfactorily addressed. 



5. Conclusion 

In this paper, we have considered two possible flow struc- 
tures in the protoplanetary disk, which correspond to two differ- 
ent prescriptions of the turbulent Maxwell-Reynolds tensor: 

(i) A viscous prescription, commonly used in the litterature, 
which typically gives rise to the so-called meridional circu- 
lation with outward flows around the midplane and inward 
flows away from it. 

(ii) A MHD prescription, which we introduce and motivate, ap- 
propriate for magnetohydrodynamical turbulence (without 
outflows) and where velocities retain the same sign over the 
bulk of the disk's thickness. 

We have compared their effects on the dynamics of small 
solids. While, in the vertical direction, the distribution falls off 
less steeply (as a gaussian) in the MHD prescription than in the 
viscous one, it is found that the net radial flux differs little be- 
tween the two. In fact, this radial flux is well approximated by a 
commonly used one-dimensional formula, which only depends 
on the vertically averaged turbulence parameter {a) p. We have 



shown that this can be generally expected for any likely prescrip- 
tion of the flow (and not only those considered here). 

It thus follows that in itself, evidence of outward transport of 
inner solar system material to the comet-forming region cannot 
be viewed as diagnostic of meridional circulation, whose physi- 
cal reality is still controversial. In fact, if the mass accretion rate 
is positive, this outward transport is rather to be attributed to tur- 
bulent diffusion. More generally, the net radial transport of early 
solar system material does not significantly constrain the flow 
structure of the gas over the disk's thickness. 



Appendix A: Exact results for vertically-integrated 
radial particle flux in the viscous and MHD 
prescription 

We list here exact results for different quantities allowing one 
to calculate the net radial flux of solids in the two prescriptions, 
assuming an equilibrium vertical distribution of the said solids. 



A.1. Viscous prescription 

{Dur)p„ ^ Dur{R,0) 



(ur)pp = ur(R, 0) ■ 



Q-(7?,0)c? d ^Jc^,\fz\^ 



D. dR \Q? 



(A.l) 
(A.2) 



r(l^2|f/)exp(-..(/^,0).->y 

<Vdrift,R>pf = Vdrift,R(^, 0) — (A.3) 



X^"exp(-5,(7?,0)ev^-/)dy 



I'p.con- - Drr- 



dlnH Ij 



dR 



with 

d\nH dlnH/dR 



d\nP dlnP(R, 0)/dR 



(z^f ^ jrf^^p{-S.(R,0)e>'-y')dy 
\h) X"exp(-5,(/?,0)e--/)dy 



(A.4) 



(A.5) 



(A.6) 



A.2. MHD prescription 



{Drr}p„ ^ Drr(R,0) 



a{R, 0) 



l-f(R,0) 



erf y^S,{R,0)\n(a,^.,Ja(R,0))) 

f(R, 0) ^ == - 

^JS,{R,Q) 

erf ( y/(S,(R,0)+l)ln{a^,Ja(R,0))) 



{ur)p= UR(R,0)f(R,0) 



y/S,iR,0)+l 



erf(V5,(«,0)ln(a„,Ja(/?,0))) 



:a.7) 



y/S,(R,0) 
_^2A (a(R,0) 



S,{Rfi)-\ 



^ffmaxc? d\na{R,Q) . 
^JnD. oR 



1 



\n{a^,Ja{R,0)) 
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B + 



a(R, 0) 



a„ 



^ 2 / aiR,0) \^'^''-°^(aiR,t 



+2A 



a{R, 0) 



S;(R,0) 



(A.8) 



dlnL 2 la(R,0) 



OR 



(A. 14) 



(Vdrift,R)p = Vdrift,R(/?,0)/(/?,0) 
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the qualitative difference between the two flows. 



2B dlnH 
+ —;= + 



d\nP 



^JS,(R,0) 

( 1 /erf ( y/S,iR, 0)ln (a^^JaiR, 0))) 

sJrMx 



^|S,(R,0) 



_2 / In iam^/ajR, 0)) ( a(R, 0) \ ^ ( "^^^ ] 

^max / 



S,(«,0) 

c 



9) 



''p.corr 



dlnf(R,0) 
dR 



(DrrY 



Pp 



-Dj^(0)S,(R,0)mO) 

1 



(dlnl. ^ d\na(R,Q) 



\ dR dR 
( erf ( y/S,iR,0)ln(a,-,,Ja(R,0))) 



S,(R,Q) 



y/S,(R,0) 



^ ^ ]n(a^/aiR,Q)) I 



a(R, 0) 



Sj(R,0) \ 



+ln 



51n2 2 laiR,0) 



amax \laiR,0)\'^^''-°^-' 2A 



a{R,G)j\ a, 



dR y/n\ a, 



S,(R,0) 



B-A 



a(R,0) 



d\nH( 1 



dR [S^iR,0) 



erf ( ^JS,(R,0)]nia,^/aiR,0))) 



_2 / In iam^/ajR, 0)7 / a(R, 0) ^ 

1 ^ \ Q^max / 



AC la(R,0) 

^Jn \ amax 



5.(R,0) \ 



(A.IO) 



with 



A = r exp {s,{R, 0) (l - ^,/) - yA dyCA.! 1) 

B = f exp (s^iR, 0) (l - \\ dy (A.12) 

C = r f.Js.(R, 0) (l - "^eA] dy (A.13) 



—lnf{R,0)=f(R,0)S,(R,0) 
dR 



Idlnl. d\naiR,0) 
\dR dR 



Vln(a„ax/a(j?,0)) / g(j?,0) \^'^"-"^^^ 
>^(5,(7?,0) + 1) \ a^ax / 

erf ( V(5z(-R,0) + l)ln(ff(/?,0)/a„ax)) 
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